Glass rheology: From mode-coupling theory to a dynamical yield criterion 
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The mode coupling theory (MCT) of glasses, while offering an incomplete description of glass 
transition physics, represents the only established route to first-principles prediction of rheological 
behaviour in nonergodic materials such as colloidal glasses. However, the constitutive equations 
derivable from MCT are somewhat intractable, hindering their practical use, and also their inter- 
pretation. Here we present a schematic (single-mode) MCT model which incorporates the tensorial 
structure of the full theory. Using it, we calculate the dynamic yield surface for a large class of 
flows. 



The 20th Century saw formidable advances in the sub- 
ject known as theoretical rheology - whose aim is to pre- 
dict or explain the nonlinear flow behavior of materi- 
als. Ideally for each class of material one wishes to gain 
a 'constitutive equation' which predicts the stress ten- 
sor at time t as a functional of the strain tensor at all 
earlier times (or vice versa). There are two broad ap- 
proaches to this task. The more traditional one focuses 
on symmetry, conservation, and invariance principles (of- 
ten of some subtlety) and then proposes empirical equa- 
tions that respect these principles p|. In the second ap- 
proach, the goal is to start from a first-principles analysis 
of molecular motion, and then by judicious (though pos- 
sibly uncontrolled) approximation arrive at a continuum- 
level constitutive model. This is clearly far more ambi- 
tious, and success has so far been restricted to relatively 
few classes of material. Perhaps the most striking suc- 
cess has been the Doi-Edwards theory for solutions and 
melts of entangled linear polymers 0, 0| (extended later 
to branched [If or breakable [f| chains). In their resting 
state, such polymers are ergodic and therefore attain the 
Boltzmann distribution: moreover their local structure is 
weakly perturbed from this, even under flow. 

Glasses at rest, in contrast, are nonergodic on exper- 
imental timescales. This poses major obstacles to the 
rheological theory of glasses, and is responsible for ag- 
ing and other phenomena that have been partially ad- 
dressed using mesoscopic models Q . The onset of arrest 
at the glass transition is, familiarly, accompanied by the 
onset of an elastic modulus. Window glass is a brittle 
solid: it deforms elastically for low stresses but shatters 
under large ones. However, some other glasses - most 
notably in colloidal suspensions (whose glass transition, 
for hard spheres, is found experimentally at about 58% 
volume fraction) are not brittle solids but show contin- 
uous yielding behavior. Although experiments suggest a 
more complex picture 0, H[ , the simplest explanation is 
that, above some yield stress, the glass melts. If a steady 
stress above the yield level is maintained, the resulting 
fluid can be expected to attain an ergodic (though non- 



Boltzmann) steady state. 

This restoration of crgodicity under steady flow offers 
one motivation for an approach to glass rheology based 
on mode-coupling theory (MCT). In particular, it miti- 
gates a well-known shortcoming of MCT which can ad- 
dress only the approach to the glass transition from the 
liquid side, and therefore cannot access aging phenom- 
ena. For systems at rest MCT predicts a true glass tran- 
sition (that is, a divergent relaxation time) rather than 
one cut off by activated processes. MCT does not ad- 
dress such processes, at least in its simplest form (for 
a discussion of extensions that do, see [1, Qlj); on the 
other hand, despite its mean-field character, MCT does 
appear to capture some aspects of dynamic heterogeneity 
fll| . When applied to colloid rheology, MCT addresses a 
set of Brownian particles advected by a flow that inter- 
act solely by conservative interparticle forces. As such it 
neglects hydrodynamic interactions whose role near the 
glass transition is unclear. (Such interactions are impli- 
cated in shear thickening, which can also arise in dense 
colloids; see e.g. fl2|.) 

The application of MCT to systems under flow re- 
quires additional approximations beyond those of qui- 
escent MCT. The most comprehensive route involves an 
'integration through transients' (ITT) formalism which 
has been detailed in a series of papers addressing in 
turn steady shear (l3| . time-dependent shear [3] and 
arbitrary homogeneous, time-dependent (incompressible) 
flows [15| ■ Below we very briefly outline the final consti- 
tutive model of (lo| . before turning to our main purpose 
in the present work. This is to present a radically sim- 
plified model, directly inspired by those results, but far 
more suitable for practical calculations of flow behavior. 
It is derived, in essence, by suppressing all wavevector 
indices on the MCT equations for coupled density fluctu- 
ations, resulting in a single-mode description. Our model 
thus stands in relation to MCT-ITT as the well known 
'schematic' models of the MCT glass transition stand in 
relation to quiescent MCT [l(|. It builds on previous 
schematic models for MCT rheology [HI, [l?} which were 
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however restricted to steady shearing, and contain no 
tensorial information about the stress response to more 
general flows. In the present work we overcome both 
of these limitations, gaining a practicable but still mi- 
croscopically motivated constitutive description for glass 
rhcology. 

In the MCT-ITT approach, the system is initially as- 
sumed to be at rest and in Boltzmann equilibrium. The 
flow is subsequently switched on; its effect via particle ad- 
vection is to create a nonequilibrium Smolochowski op- 
erator. An exact nonequilibrium Grccn-Kubo relation 
allows the stress tensor at any later time to then be ex- 
pressed as a time integral of the expectation value of a 
certain operator evaluated within the equilibrium ensem- 
ble [13, E3- Applying MCT-type factorizations to this 
expectation value gives a constitutive description in three 
parts. The first relates the deviatoric stress to strain his- 
tory via a time integral whose (wavevector dependent) 
kernel involves derivatives (with respect to wavenumbcr 
k) of the equilibrium static structure factor Sk and the 
square of a normalized transient density auto-correlator: 
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This expression is found from results reported in [15| 
(suppressing an isotropic pressure term) and is a con- 
sequence of approximating the relaxation of stress fluc- 
tuations by that of the density fluctuations responsible 
for slow stuctural relaxation in dense systems. <I>k(£, t') 
is a transient density correlator which, by virtue of the 
MCT-ITT approach, is calculated using the equilibrium 
distribution. The correlator measures the overlap of a 
density fluctuation at wavevector k and time £ with one 
at earlier time £' whose wavevector k(£, £') evolves due to 
flow-induced advection to become k at £. For the pur- 
poses of calculating this correlator, the system is taken to 
be in equilibrium at £', with the flow acting thereafter. 
Because averaging in MCT-ITT can be done with the 
equilibrium distribution function, only the equilibrium 
structure factor Sk appears in Eq.l. In Eq[TJ B is the 
Finger tensor: 



B(t,t') = E(t,t') ■ E T {t 1 t'), 



(2) 



where E(t,t') is the deformation tensor applied between 
a past time £' and the present time £. (This obeys 
dhxE/dt = n(t), where k is the strain-rate tensor.) 
The reverse-advected wavevector is given as k(£, £') = 
k • E(t,t'). (Note that, where these differ, we use the 
notation of 0] rather than Q. Note also that we follow 
the notation and formalism of [l5| which is modified from 
that of [H Q [13, [3 as explained in detail in 0.) 

The second part of the MCT-ITT description deter- 
mines the time evolution of the correlators appearing in 



Eq[T]in terms of a three-time memory kernel: 
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where the ovcrdots denote partial differentiation with re- 
spect to the first time argument. Here the 'initial decay 
rate' obeys r q (£, £ ) = D q 2 (t, t )/S^ t ,t ) with D a bare 
diffusivity, and q(t,t') = q ■ Finally, MCT- 

ITT approximates the memory kernel m q (£, £', £o) by the 
factorized expression 



m q (£,£',£ )=- 



(4) 



qv™, 167r3j r~ g2(£',£ )g2(£,£ ) 
X Vq kp (£' , t ) ^qkp (t, *o)*£(f ,t ) (t, t')$ P (t< , t0 ) (t, 0, 
where p = q — k, and the vertex function obeys 
Vqkp(i,to) = q(*,to) • (k(£,£o)cfc (Mo) +p(£,£o)cp( f ,t )), (5) 

with Cfc = 1 — 1 / Sk- Again, like in Eq[TJ this assumes that 
fluctuating stresses decay mainly by density fluctuations. 
In EqsflHSl MCT-ITT offers a closed constitutive model 
for interacting Brownian particles near the glass transi- 
tion, requiring only Sk and k(£) as input to calculating 
the stress. 

So far, however, these equations have proved amenable 
only to approximate solution even for the simplest flows, 
such as steady and stepwise deformations [13 . 0, Ell- 
in three dimensions, at least, such approximations have 
all involved replacing the complicated angular depen- 
dendence of correlators in k-space by an isotropic de- 
pendence, creating the so-called isotropically sheared 
hard sphere model or 'ISHSM' [l3]- (Avoidance of 
this has however very recently proved possible for two- 
dimensional systems [20j| .) Our aim here is to create a 
much more tractable, simplified description that nonethe- 
less captures both the tensorial character and the basic 
time-dependence of the full MCT-ITT equations. In do- 
ing this we pay due attention to invariance properties 
which are respected by the microscopic theory but which 
could get lost in approximation. Of these, the least ob- 
vious is the principle of material objectivity, which re- 
quires invariance of the constitutive model under arbi- 
trary time-dependent rotations. This stems from neglect 
of inertia (implicit in our description based on Brownian 
particles) and imposes strong constraints on the tensorial 
structure of rheological models [3] ■ 



A. Derivation of schematic model 

Our starting point is EqQ]for the stress tensor; to cre- 
ate a single-mode description we discard all wavevector- 
dependent information. A first step is to discard all an- 
gular information. To do so we assume that both the 
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correlator $ and the advected wavenumber k(t, t') are 
isotropic in k, so the integrand in Eq[T] becomes a prod- 
uct of an anisotropic and an isotropic term (each en- 
closed in square brackets). To see that this factorization 
respects the principle of material objectivity, note that 
the advected wavenumber k(t, t') is not altered by a time 
dependent rotation Q(t), while the Finger tensor in the 
rotated frame is given by B{t,t') = Q(t)B(t,t')Q T (t). 
Substitution into Eq[T] and changing integration variable 
to k' = k • Q{t) yields &(t) = Q(t) a(t)Q T (t) as desired. 

With the isotropic assumption we can perform the an- 
gular integral in EqfT] explicitly, to obtain 
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Having dealt explicitly with the tensorial structure, we 
may now safely discard all wavevector dependence in the 
expression for the generalized modulus G(t, t'). Thus we 
obtain a fully schematic model in which Eq[5] holds with 



G(t,t') = v a <£ 2 (t,t') 



(7) 



where v a = G(t, t) is a parameter measuring the strength 
of stress fluctuations. When analyzing flow curves of hard 
sphere-like colloids values of the order lOOfceT / d 3 are 
typically obtained, where d is the sphere diameter [2lj . 

The full microscopic equation of motion for the cor- 
relators is given by Eqs[3][5] Discarding all wavevector 
indices in Eq[3] gives the schematic representation 
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There is now only one initial decay rate, T, which sets 
the microscopic time scale and may thus be set equal to 
unity. 

The microscopic form for the memory function (Eqf5|) 
depends quadratically on the correlators and has ver- 
tices (Eqj5]) that are a function of the accumulated strain 
over different time intervals. This strain dependence of 
the memory kernel is a central novel aspect of our mi- 
croscopic constitutive equation (Eqs[T]-[5]) the qualitative 
aspects of which we aim to reproduce in a simplified 
schematic model. In the absence of flow, our memory 
function reverts to that of standard (quiescent) MCT, 
with a time-independent coupling to a product of two 
stationary correlators. The arrest transition of this stan- 
dard MCT is captured by a minimal schematic model 
(replacing in effect the k integral by a single wavevec- 
tor) known as the F\ 2 model [l(|, with both a linear 
and a quadratic coupling. It is well established that the 
schematic i 7 " 12 model captures quantitatively many uni- 
versal aspects of the glass transition contained in the full 
MCT equations. To extend the F\ 2 model to address 



rhcology we now make the ansatz 

m(t,t , ,t ) = h 1 (t,t )h 2 (t,t')h s (t , ,t ) x 

whose form is inspired directly by Eq[H So long as the 
product of the three h factors becomes constant in the 
absence of flow, the standard F\ 2 model is recovered in 
that limit. 

Within the microscopic description, the wavelengths 
of density fluctuations are on average reduced by defor- 
mation, allowing diffusion and intcrparticle forces to re- 
lax these fluctuations more rapidly, so that strain causes 
a progressive loss of memory (this mechanism is repre- 
sented by the time dependence of the vertex functions in 
Eq.4). The hi should thus be decaying functions of the 
strain accumulated between their two time arguments. 
However, there is some redundancy between these fac- 
tors and we have found that if hi and h 2 both decay, /13 
can be set to unity without losing any obvious physical 
content of the model. In addition we set hi = hi = h 
so that the model contains only one unknown function 
h(t\, t 2 ); choosing this as detailed below gives good qual- 
itative agreement with the full theory for all flows so far 
investigated. By these considerations we arrive at the 
following expression for the schematic memory function: 

m{t,t',t ) = h{t,t )h(t,t') [ui$(M') + v 2 $> 2 {t,t')]. (9) 

Our final task is to decide the form of h(ti,t 2 ). This is 
shorthand for h(E(ti,t 2 )), and encodes the loss of mem- 
ory caused by strain. Earlier work on steady shear shows 
that the choice ft. = 1/(1 + (7/7c) 2 ) (with 7 = (t 2 — £1)7 
the accumulated shear strain) gives results close to the 
full MCT-ITT theory for this case (T3l ■ We can be guided 
by this finding, but need to generalize it to non-shear 
flows, such as elongational deformation, which can be ex- 
pected to differ quantitatively in their effects on memory 
loss. Notably, this variation is strictly limited in form: 
for incompressible systems, the principle of material ob- 
jectivity ensures that the deformation tensor E can enter 
only through the invariants I\ = Tr£? and I 2 = (Trf? -1 ) 
||. We are thus led to the following choice 



h(t,t ) 



7c 2 + [uh{t,t ) + (1 - v)I 2 (t,to) - 3] 



(10) 



where we have introduced both a mixing parameter v 
(0 < v < 1) and a crossover strain parameter j c . (The 
latter sets the scale for the recoverable elastic strain; with 
this as a fit parameter, schematic models accurately ac- 
count for structural distortions in steady state shear [2l|.) 
Note that I\ = I 2 for both shear and planar extensional 
flows whereas I\ 7^ I 2 for uniaxial elongation. 

EqsISl flTjl fully specify our new schematic model. In a 
further (conventional) simplification, the parameters V\ 
and v 2 are now replaced by v 2 = 2 and v\ — 2(y/2 — 1) + 
e/(v2 — 1). The 'separation parameter' e is analogous 
to the density in a microscopic system, with negative 
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values corresponding to fluid states and positive values 
corresponding to glass states. We have thus obtained a 
closed constitutive model with four adjustable parame- 
ters (v ai r, 7 C and v) and one control parameter (e). This 
compares with the full MCT-ITT equations in which the 
static structure factor S q , density and bare diffusivity 
control all measurable quantities. Accordingly one can 
view tv, 7 C and v as fixed in some unspecified way by 
S q . The latter in turn depends on the thermodynamic 
control parameters and on the interaction potential be- 
tween particles. When applied to experimental data, all 
five parameters of the model are used to fit the data and, 
as expected, exhibit smooth variation with the external 
parameters pH ]. 



B. Physical content of the model 

The constitutive model laid out above describes a fluid 
whose instantaneous elastic response is linear in the Fin- 
ger tensor B, as manifest in Eq[6] The Finger tensor is a 
natural frame-invariant extension of the linear strain ten- 
sor, and a similar dependence of stress on strain arises, 
for instance, in models that involve networks of Hookcan 
springs. Indeed, were G(t,t') in EqJS] to be replaced 
by a strain-independent, time-translation- invariant ker- 
nel G(t—t'), this would recover the Lodge equation [l|,[H, 
a standard rheological model. (The relationship to such 
models is explored further in the Supplementary Informa- 
tion, SI.) The additional nonlinearity in our model arises 
almost entirely from the strain-induced erasure of mem- 
ory. This differs somewhat from EqQ] in which there is 
some additional nonlinear elasticity, even setting $ = 1; 
but this does not survive the schematic-model relaxation 
approximation that gives Eq[71 

Thus the kernel G(t,t') in EqJ7]is strongly sensitive to 
the strain-induced decay of the correlators, as found in 
turn via Eq|3] This decay sets in at strain increments of 
order 7 C and is absent in the linear response regime for 
which G(t - t') = v a $ e (t - t') 2 with <f> e the quiescent 
state correlator. For e < linear response describes a 
viscoelastic fluid; for e > 0, a viscoelastic solid. In the 
latter (glass) regime, a sustained strain rate causes the 
memory (and therefore the correlators) to fall to zero 
at long times, resulting in fluidization. When such a 
flow ceases, the system solidifies again. Because of the 
memory-function structure, the system's current mate- 
rial properties (for example, the relaxation time that gov- 
erns response to an additional strain increment) are in- 
fluenced by past deformations; but this memory is itself 
erased by large strains. Thus, for instance, if a steady 
flow is suddenly switched off, it takes some time for the 
memory kernel to rebuild, the correlators to slow down, 
and the solid-like properties of the glass to be reinstated. 

Note that the model captures memory erasure by sud- 
den (step-strain) events as well as by sustained deforma- 
tion rates. Here there are considerable subtleties to the 
model. For instance m(t,t',to) is diminished by a step 
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FIG. 1: The steady state shear viscosity n (full lines) and first 
normal stress coefficient = (a xx — Cyy)/^ 2 as a function 
of Peclet number Pe for three values of e. At low shear rates 
both n and ^i saturate to a plateau value for fluid states 
but diverge in the glass. (Parameters: r = 1,7 C = 0.1, v a = 
100,^ = 1.) 



strain that occurs at any time between t and t, as makes 
sense given that the time integral in Eq[5] spans this in- 
terval. However, if this step strain is later reversed, only 
pairs of steps that bracket the intermediate time t' in that 
integral contribute to a reduction of memory. Such fea- 
tures are sharpened by taking the schematic limit, but so 
far as we can tell they do have proper anticedents within 
the full model, Eqs[3][5] 



I. RESULTS 

As mentioned above, in steady shear our tensorial 
schematic model gives for the flow curve o- xy (f)) results 
similar to those of a model (there called the F± 2 model) 
published previously [13, HH, H3| • The latter model how- 
ever says nothing about other elements of the stress ten- 
sor; nor nonsteady flows (23|; nor non-shear flows. Below 
we give sample results from the new model that illumi- 
nate each of these issues in turn. Unless otherwise stated 
we select parameters r — 1, 7 C = 0.1, v a — 100, v = 1. 
The first of these choices, T = 1, sets the time unit; our 
choice directly equates strain rate to the Peclet number, 
Pe. (Continuing to assert this relation while choosing 
a different T would allow flow curves to be scaled hori- 
zontally, possibly improving the fit to experiment.) The 
values chosen for 7 C and v a reflect the values seen exper- 
imentally for these two quantities (~ 0.1, 100, with the 
latter in units of fcsT/df) in hard sphere systems close 
to the glass transition [2J|. Note that smaller or larger 
values might be appropriate in systems with non hard- 
sphere interactions, such as short-range bonding forces, 
which, in the full model, enter through changes in Sk- 
Turning to the mixing parameter v, our results for sim- 
ple shear and planar extcnsional flows are independent 
of this parameter and for other flows display only a very 



weak v dependence. For simplicity we therefore choose 
v = 1. The remaining variables are s, the distance to the 
glass transition, and the flow history. 

A. Normal stresses under steady shear 

The tensorial schematic model predicts the following 
diagonal stress tensor elements under steady shear flow 

/>oo 

a xx = 2fv a / dt't'$ 2 (t') , (11) 
Jo 

Oyy = a zz = 0, 

where the correlator only depends on a single time ar- 
gument due to restoration of stationarity. Defining the 
first normal stress coefficient as = (a xx — <Ty y )/j 2 = 
iVi/7 2 , we find 

poo 

*i = 2w a / dt't'$ 2 (t'). (12) 
Jo 

Here the correlator is the one calculated with the flow 
present, which decays to zero at long times. Thus is 
finite and positive, meaning that in a steady shear ex- 
periment the bounding plates tend to get pushed apart. 
A positive ^1 has been found recently in a full numeri- 
cal solution of the MCT-ITT equations in two spatial di- 
mensions [1^ , and the same is found for the full model in 
three dimensions if one invokes the isotropized (ISHSM) 
approximation. On the other hand, the second normal 
stress difference N2 = cr yy —a zz is, from Eq llll identically 
zero, and this is also true of ISHSM. 

FigQ] shows steady-state flow curves for the viscosity 
and the first normal stress coefficient in simple shear from 
the present model. In the fluid phase, e < 0, both follow 
a power law as a function of the separation parameter, 
with 77(7 0) ~ (-eT 2 - 34 and ^(7 -f 0) ~ (-e)" 4 - 68 . 
For e > (within the glass phase) both the viscosity 
and the first normal stress coefficient remain divergent 
as 7 — > 0, consistent with the presence of a dynamic 
yield stress - that is, a finite limiting stress tensor at 
low shear rates. (As shown in [l7j], the full MCT-ITT 
predicts a nonzero dynamic yield stress throughout the 
glass phase; so does our schematic model as detailed be- 
low.) As mentioned in the derivation of EqJHl the viscos- 
ity curve almost coincides with the model of p7| . 
which in turn closely approximates the full ISHSM pre- 
dictions, and which themselves give a good account of 
experimental data [U, [24[. Note that, with the choice 
of units conventionally adopted in MCT-ITT, the Peclet 
number in Fig[T] and subsequent figures coincides with 
the strain rate, Pe= 7c? 2 /Dq with Dq the (bare) diffusiv- 
ity. 

B. Step shear strain 

In the case of a mathematically idealized step shear 
strain, j(t) = 7<5(t), EqJT] simplifies considerably (To 
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FIG. 2: The long time shear stress following a step strain as a 
function of the step amplitude 7 for three glassy statepoints. 
At small values of 7 the expected linear response behavior is 
recovered. For larger 7 values a regime of sublinear increase 
indicates the onset of plastic flow. (Parameters: T — l,-y c = 
0.1,1V = 100,i/= 1.) 



attain this limiting case, one requires a strain ramp of 
duration t -C which may not be achievable experi- 
mentally. We leave the study of finite ramp rates to fu- 
ture work.) The subsequent shear stress is given simply 
by 

O xy {t) = V al <$> 2 {t), (13) 

where the correlator <f>(i) satisfies the linear equation 

+ r($(i)+ / d£ m{t-t')$(t')) = 0,(14) 
Jo 

Here ^eq(i) is the equilibrium correlator obtained by 
solving Eq[8] in the absence of flow. The first normal 
stress difference N± = a xx — a yy for times following the 
step is given by 

iVi(i) =v al 2 ® 2 (t). (15) 

This satisfies the Lodge-Meissner relationship 
Ni(t)/a xy (t) = 7 which holds for all constitutive 
models in the (misleadingly named) 'simple fluids' class, 
of which ours is a member [l[. 

The shear stress after step strain relaxes monotoni- 
cally to zero in a fluid but has a nonzero asymptote 
ct(oo) > in the glass. Numerical calculations of ct(oo) 
using the ISHSM [lj] showed this to be a nonmonotonic 
function of strain amplitude. This is a consequence of the 
strain-erasure of memory which can over-compensate the 
(essentially linear) dependence of the initial stress level 
on the strain, causing curves for different strain rates to 
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FIG. 3: The planar extensional Trouton ratio for two fluid 
statepoints (s < 0) and one glassy state (e > 0) calculated 
using the schematic model (full lines). For comparison we 
show also the microscopic results from [l5| (points connected 
by lines) where the full model was solved in the isotropic 
(ISHSM) approximation. As the glass is approached the lin- 
ear response regime vanishes, leading to a nontrivial Trouton 
ratio. The Peclet number employed in the microscopic cal- 
culations is denoted Pe m i c . (Parameters: T — 10~ 3 ,7 C = 
0.1, v a = 100,i/= 1.) 

cross. The physics of this effect is captured within our 
schematic model; Figj2]shows ct(oo) as a function of 7 for 
various e and this is quite similar to Fig.l of There 
are deviations at large strain, caused by our choice of 
Eq[l0]to describe the strain effect. In fact, our schematic 
model may be more physical here than ISHSM which 
is found numerically to predict negative cr(oo) at large 
enough strains [ljj]. (It is not yet clear whether the lat- 
ter is an artefact of the isotropization used or signifies a 
deeper problem with MCT-ITT itself.) 

C. Steady planar and uniaxial elongation 

For these flows, as in steady shear, our schematic 
model produces qualitatively similar flow curves to 
ISHSM. A sensitive test of the tensorial aspects of the 
model is to plot the "Trouton ratio", (a xx — cr yy )/a X y 
as a function of Pe. Fig[3] shows a direct compari- 
son for the case of a planar extensional flow (fi^ = 
j(5 X iS m j — Syidyj)). For a fluid phase this must asymp- 
tote to 4 at small strain rates, but in the glass a nontriv- 
ial (and smaller) value is possible. Within the schematic 
model, this nontrivial value is numerically indistinguish- 
able from 2, while the ISHSM gives a rather smaller num- 
ber; for uniaxial extension the schematic model gives a 
noninteger value (close to 1.7), this time slightly above 
ISHSM (~ 1.6). A quantitative comparison of this kind 
requires a parameter-matching exercise to be undertaken. 
For this purpose, choosing e = 3(0— <j) c ) roughly matches 
the correlator's final relaxation time to a typical one in 



the ISHSM at volume fraction <fi [l7|. To also match the 
rate of change of the yield stress with e to its ISHSM 
counterpart would require additional parameters in the 
schematic model. Qualitatively, this difference can be ab- 
sorbed into an effective Pe number (equivalently, T ^ 1) 
as in FigJ3] Once this is done, the Trouton curves for 
schematic and ISHSM results are qualitatively very sim- 
ilar. 



D. Yield surface 

One striking feature of MCT-ITT is that it allows in 
principle an unambigous first-principles determination of 
a dynamic yield stress for glasses. For any given flow 
geometry (e.g. simple shear, planar elongation, or uni- 
axial elongation) a yield stress tensor is defined as the 
limiting stress obtained as the relevant flow rate tends to 
zero from above. One expects that the manifold of these 
yield stresses divides flowing states at large stress from 
non-flowing states at small stress and that this manifold 
forms a closed surface in some appropriately constructed 
space. 

We emphasize that the dynamic yield manifold may 
differ from that for static yield. The latter is usually 
defined as the step stress that must be exceeded to re- 
sult in steady flow; this may depend on prior sample 
history whereas the dynamic yield manifold (being de- 
fined as the limiting stress within a series of ergodic, 
fluidized steady states) does not. The relation between 
the two yield surfaces directly mirrors the familiar one 
between dynamic and static friction, and although they 
address different physical situations, each merits care- 
ful study. Here we present results only for the dynamic 
case which, because of the simpler flow history, is easier 
to compute. Crucially, however, our schematic model is 
tractable enough to allow prediction of the long-time re- 
sponse to step strains (or step stresses) of arbitrary char- 
acter and amplitude (Fig[2|). Combining these, we may 
compile a conventional static yield manifold, and also ad- 
dress many other aspects of the flow-history-dependent 
yield response. We defer to future work an exploration 
of these problems, but believe that our model offers a 
promising new semiquantitative route, inspired by first 
principles statistical mechanics, into the study of diverse 
yield phenomena in amorphous solid materials, including 
static yield. 

To calculate the dynamic yield manifold, we note first 
that for any particular flow geometry, the limiting stress 
tensor at vanishing flow rate can be diagonalized to give 
principal stresses s\,S2,s^. Eliminating the (arbitrary) 
isotropic pressure for an incompressible system, we can 
then plot a point in the two dimensional (s\— S2 , S2—S3) 
plane. Repeating this procedure for all possible flows 
leads to a locus of points that define the 'yield stress 
surface' in this two dimensional space. We call this 'rep- 
resentation 1' of the yield manifold. In an alternative 
representation, all three stress variables are retained to 
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give a surface in a three dimensional space; but since hy- 
drostatic pressure is irrelevant (due to incompressibility) 
this surface must be translationally invariant along that 
axis. Thus looking down the hydrostatic axis again gives 
a closed curve (representation 2). 

Using symmetry and other arguments, it is possible 
to map out an entire yield surface by considering only a 
one-parameter family of flow geometries that interpolate 
between planar and uniaxial flows. (The explicit con- 
struction is given as Supplementary Information). The 
case of simple shear is not among this family of flows, but 
its stress is among the family of stresses spanned by the 
family, and shear flows in fact lie on the same yield sur- 
face to numerical accuracy. The resulting yield surface 
is shown, in representation 2, in FigSl This looks sug- 
gestively circular, but on close numerical study is found 
to deviate disccrnibly from a circle at the percent level, 
with the maximal deviation at the points of uniaxial ex- 
tension. Our results for the yield surface are independent 
of r and are insensitive to variations in v. Increases in 
e, v a and 7 C simply lead to a scaling of the surface, as 
the yield stress increases, e.g., on moving deeper into the 
glass. 

Although ours is not a first-principles approach, but 
merely a schematic model inspired by one, we think the 
calculation of this yield surface remarkable for two rea- 
sons. First, the yield surface in representation 2 is almost 
perfectly circular. Perfect circularity (corresponding to 
an ellipse in representation 1) is the content of 'von Mises 
law', an empirical relation which has been used for nearly 
a century to approximate the yield surface of a wide range 
of materials [25[ . The yield surfaces in question are static 
ones, whereas we calculate the dynamic: we therefore 
make no claim to have explained von Mises law, but do 
find the correspondence intriguing. It will be very inter- 
esting to see whether the static yield manifold computed 
from our model has a similar degree of circularity. In any 
case, since anisotropy must cause the yield condition to 
depend on material axes as well as principal stress axes, 
the starting point for any statistical-mechanical 'deriva- 
tion' of von Mises must presumably involve an isotropic 
(and therefore amorphous) solid - that is, a glass. Our 
work represents a promising starting point for such an 
endeavor. 

Secondly, the yield surface is not quite circular, de- 
spite all the simplifications that have been made in de- 
riving the schematic model. Indeed it can be shown for 
the schematic model that the circular yield manifold is 
recovered only under conditions where the first normal 
stress difference, at the point of yield under simple shear, 
is negligible (see Supplementary Information for details). 
This limit is approached whenever yield strains are small 
and might form the basis of a systematic expansion with 
a von Miscs-like circle as the zeroth order contribution. 
As found in experiments [24j the yield strain for hard 
sphere colloids is indeed modest. 
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FIG. 4: The dynamical yield surface for a glass with e = 10 
as viewed along the hydrostatic axis sj = S2 = S3 in the space 
of principle stress differences. The red points correspond to 
planar extensional flow and the blue points to uniaxial ex- 
tensional flow. The surface is not perfectly circular and ex- 
hibits maximal deviation from the von Mises circle at points 
of pure uniaxial extension. (Simple shear flows lie too close 
to planar extension to be shown separately.) (Parameters: 

r = 1, 7 C = 0.1, w CT = 100, v = 1.) 



E. Conclusion 

We have presented a schematic model for the consti- 
tutive rheology of glasses, suppressing all wavevector in- 
dices but retaining much of the tensorial content of the 
full MCT-ITT approach 0, [3, EH . The resulting tenso- 
rial structure of the schematic model satisfies applicable 
invariance laws for the nonlinear flow of materials with- 
out inertia [2]. Our schematic model performs well at 
capturing qualitatively the behaviour of the full theory. 
More precisely, it gives results that are similar to those 
found by the (wavevector-dependent but still simplified) 
ISHSM approximation, in all flows for which results for 
the latter have so far been obtained. By construction, 
however, the schematic model is much easier to imple- 
ment numerically across a wider range of flow geometries 
and histories. As a concrete example of this, we have 
calculated the full dynamic yield stress manifold by ad- 
dressing a family of steady flows that interpolate between 
planar and uniaxial elongation. The resulting yield sur- 
face is very similar, but not quite identical, to the empiri- 
cal form of von Mises that has been widely used to model 
(static, not dynamic) yield and plasticity in solids. It re- 
mains to be seen whether the static yield manifold for 
our schematic model is also of this form. More generally, 
our work offers promise for a better understanding of the 
physics of plasticity based on statistical mechanical prin- 
ciples applied to amorphous, isotropic solids - specifically 
glasses - and is a step towards the rational prediction of 
the general nonlinear rheology of this important group of 
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materials. 
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Appendix 

Relationship to Lodge and Maxwell equations 

The present schematic model is closely related to the 
Lodge equation Q 



*{t) 



dt' B(t,t')^—- , (16) 

-r^i A 



where A is the relaxation time of the system, G is the 
instantaneous shear modulus and B(t,t') is the Finger 
tensor measuring the relative strain. The Lodge equa- 
tion is valid in both the linear and nonlinear regimes and 
is a standard model in continuum rheological modelling, 
derivable from a number of simple molecular models, e.g. 
the dumbell model for dilute polymer solutions Q ■ Inte- 
grating Eq |16l bv parts yields the alternative form 

<r(t) = fjt> (-A B (t,f)) Ge-t-W (17) 

which has a structure close to that of our schematic Eq.6. 
In this spirit we integrate Eq.6 by parts to obtain a 
Lodge-type expression 



<r{t) = v a J t dt' B(t,t')^$ 2 (t,t') 



(18) 



The present theory thus goes considerably beyond the 
standard Lodge equation by incorporating memory which 
is both nonexponcntial and a function of two-time argu- 
ments, reflecting the lack of time-translational invariancc 
under general flow. Replacing the correlator with a sim- 
ple exponential trivially recovers the Lodge equation. 

The differential form of the Lodge equation is the 
upper-convected Maxwell equation [3j , which is a simple 
nonlinear generalization of the familiar Maxwell model of 
viscoelasticity [26J. It is thus of interest to consider the 
differential form of Eq.6. Differentiation of EqfT8] yields 
an integro-differential equation for the stress 



Do- 



d 



dt' B(t,t')[—$%t')) = 21^,(19) 



where T is the initial decay rate appearing in Eq.8, 6 
is the unit tensor and we have introduced the upper- 
convected derivative [23 



(20) 



If we neglect particle interactions and therefore memory 
effects in Eq.8 the correlator exhibits pure exponential 
decay, = cxp(-r(t - t')), and EqfT9l reduces to 

the upper-convected Maxwell equation Q 



A- 



Dcr 

~Dt 



(21) 



where A = l/2r. The fact that Eqfl9l cannot generally 
be reduced to simple differential form reflects the non- 
exponential character of the correlator. We note that the 
appearance of the upper-convected derivative in Eq ll9l is 
due to the tensorial dependence on B(t,t') (a prefactor 
B -1 (t, i') would yield a lower-convected derivative). 

The Lodge equation is a member of the K-BKZ class 
of equations Q , for which 



a(t) 



dt' 



dit 
dh 



B{t,t') 



Oil 

dh 



B~ l {t, t') ] . (22) 



with u(I\, I2, t — t') a function whose time integral gives 
the strain energy. This class of equations first arose in the 
context of rubber elasticity theory and depending on the 
choice of Kernel u admits a great variety of behaviour. 
Even allowing for breakdown of time translation invari- 
ance, so that u = u{I\, I2, t, t'), we find that our model 
is not of this form (except for the limiting case of v = 1). 
However, modulo the same generalization of the time ar- 
guments, our model does belong to an extension of the 
Rivlin-Sawyers class of equations [27| in which the coef- 
ficients of B and B~ x in Eq[22] need not satisfy an in- 
tegral relationship. This more general class of equations 
are capable of predicting a richer variety of rheological 
phenomena than that of the K-BKZ class. We note that 
the generalization to two time arguments in the strain 
energy kernel u has nontrivial consequences and can give 
rise to a stress response significantly different from that 
of the standard K-BKZ theory (even in the case v = 1). 



The von Mises criterion 

The von Mises criterion is a well established empirical 
rule for determining the onset of plastic yield under an 
applied stress 25| . The criterion is typically expressed in 
terms of the principal stresses 



- s 2 f + (s 2 S3 ) 2 + ( Sl - , S3 ) 2 ) = «) 2 , (23) 



where a v ss is the shear stress at yield under simple shear 
deformation. Eq |23l is based on the physical assumption 
that yielding occurs when the distortion strain energy 
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(defined as the left hand side of Eqj23] divided by 2G, 
where G is the shear modulus) exceeds a critical value. 
It is useful to interpret Eql23l geometrically in the space 
of principal stresses, where it describes a surface separat- 
ing elastically deformed states from states of plastic flow. 
Eql23l defines a circular cylinder with axis along the line 
Si = S2 = S3 and radius \[2a v as . The invariance along 
the line s\ = s 2 — S3 is a geometrical reflection of the 
fact that the yield condition is independent of hydrostatic 
pressure. 



Flow parameterization 

For calculation of the dynamic yield stress surface we 
employ a one-parameter family of flow geometries which 
interpolate between planar and uniaxial flow. Consider 
the following parameterized Finger tensor 



Defining the integral 



B = 




(24) 



By varying the parameter A from to 1 the flow changes 
continuously from uniaxial to planar elongation. Using 
our schematic model to calculate the dynamic yield stress 
for each value of A in this range enables one twelfth of 
the yield surface to be completed. General symmetry re- 
quirements then enable the yiel d surface to be completed 
without further calculation [25[. The flows characterized 
by Eql2H for < A < 1 are thus physically distinct. 
Other cases, such as uni-biaxial flow (A < 0), are de- 
generate in the sense that they do not provide additional 
information regarding the yield surface. 

Substitution of EqJ2H into Eq.6 yields the principal 
stresses 



si = v a j j dt' (2e 2 ^') $ 2 (0, 

s 2 = wj^dt 1 (-(1 + A)e-( 1+A ^')<i> 2 (0, 

T j J dt' (-(1 - A) e-t 1 - 4 )^') $ 2 (t'), (25) 



S3 = v Q 

where the correlators are those calculated using Eqs.8- 
10 with the flow given by EqJMl Connection to the von 
Miscs criterion can be made by expanding the exponen- 
tials in Eq |25l to first order in jt. This yields a moment 
expansion 

/>oo />oo 

—=27/ dt'<S> 2 (t') +4 7 2 / dt' Y$ 2 (i'), (26) 

v <r Jo Jo 

pOO pOO 

— =-(A+l)j dt' <P 2 (t') + {A + l) 2 j 2 dt't'$ 2 (t'), 
v* Jo Jo 



^={A-l)j dt'$ 2 (t') + (A-l) 2 j 2 dt't'&U'). 

v <y Jo Jo 



Jn = V a j 



n+1 



dt' (t') n $ 2 (t'), 



(27) 



and substituting Eqs[26]into the left hand side of Eql23l 
yields 



g((si-s 2 ) 2 + (s 2 -s 3 ) 2 + (s 1 -s 3 ) 2 ) = 

(3 + A 2 )J 2 + i(3 + A 2 ) 2 J 2 + 6(1 - A 2 ) J x J 2 (28) 
o 



Using Eqs.8 and 9 it can be shown that in the limit of 
small flow rate 7^0 



VsTa 2 " 'Jo 



7 / dt'$ 2 (t>) 



V3 + A 2 "' 



A -3T^^i' it ' t ' 4i ^ <t ' ) = PT^) (2i 



where Nf is the first normal stress difference at yield 
under steady simple shear. Substitution of Eq l29l into 
Eql26l thus yields 



- ((si - s 2 ) 2 + (s 2 - S3) 2 + (si - S3) 2 ) = 



12 



(Nf) 



y\2 



3(1- A 2 ) 
(3 + A 2 f/i 



Nfal. (30) 



This expression is the dynamic analogue of the von 
Mises criterion for static yield. Neglecting the first nor- 
mal stress difference reduces EqJ3Dl to the von Mises 
form. Numerical calculations using the schematic model 
show that Nf is approximately two orders of magnitude 
smaller than a^ s for a given value of e. The third term on 
the right hand side of EqJ3D] thus constitutes a small cor- 
rection to circularity. The dependence of the third term 
on the parameter A lends non-trivial structure to the 
yield surface. Continuing the moment expansion (Eq|26p 
to higher orders in jt generates further small corrections 
and may thus be viewed as a systematic expansion about 
a circular dynamic yield surface. The numerical results 
presented in the main text were performed using Eg 1251 
without further approximation. 
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